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We investigate the effect of spin-orbit and spin-spin couplings on the estimation of parameters 
for inspiralling compact binaries of massive black holes, and for neutron stars inspiralling into 
intermediate-mass black holes, using hypothetical data from the proposed Laser Interferometer 
Space Antenna (LISA). We work both in Einstein's theory and in alternative theories of gravity of 
the scalar-tensor and massive-graviton types. We restrict the analysis to non-precessing spinning 
binaries, i.e. to cases where the spins are aligned normal to the orbital plane. We find that 

■ the accuracy with which intrinsic binary parameters such as chirp mass and reduced mass can be 
estimated within general relativity is degraded by between one and two orders of magnitude. We 
find that the bound on the coupling parameter ljbd of scalar-tensor gravity is significantly reduced 
by the presence of spin couplings, while the reduction in the graviton-mass bound is milder. Using 

' fast Monte-Carlo simulations of 10 4 binaries, we show that inclusion of spin terms in massive black- 

On, hole binaries has little effect on the angular resolution or on distance determination accuracy. For 

■ stellar mass inspirals into intermediate-mass black holes, the angular resolution and the distance are 
determined only poorly, in all cases considered. We also show that, if LISA's low-frequency noise 
sensitivity can be extrapolated from 10 -4 Hz to as low as 10 -5 Hz, the accuracy of determining 
both extrinsic parameters (distance, sky location) and intrinsic parameters (chirp mass, reduced 
mass) of massive binaries may be greatly improved. 
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I. INTRODUCTION 



o 

The Laser Interferometer Space Antenna (LISA) is being designed to detect gravitational-wave (GW) signals in 
the frequency band between Hz and 10 _1 Hz 1]. Ground-based interferometers such as LIGO, GEO, VIRGO 
and TAMA are sensitive in a higher frequency band, between 10 Hz and 10 3 Hz. Operating at these low frequencies, 
LISA can detect, among other sources, inspirals and mergers of massive black holes (MBH) with masses in the range 
10 4 -10 7 Mq. Another promising source is the inspiral and capture of stellar-mass compact objects - such as neutron 
. £h I stars (NS) or black holes (BH) - by intermediate-mass black holes with masses in the range 10 2 -10 5 M©. 

Gravitational radiation reaction drives the inspiral of these binaries. The amplitude and phase of the gravitational- 
wave signal carry information about binary parameters, such as masses and spins, and about the location and distance 
of the binary. They may also be different in different theories of gravity. Therefore LISA can provide important 
astrophysical information, yield interesting tests of fundamental physics, and place bounds on alternative theories 
of gravity. In this paper, we consider, along with standard general relativity, theories of the scalar-tensor type 
(the simplest exemplar being that of Brans and Dicke) and theories with an effective mass in the propagation of 
gravitational waves (which we call massive graviton theories, for short). In scalar-tensor theories the phasing evolution 
is modified predominantly by the presence of dipole gravitational radiation reaction in the orbital evolution (in general 
relativity the lowest radiative multipole moment is the quadrupole). In massive graviton theories the gravitational 
wave propagation speed depends on wavelength: this generates a distortion in the time of arrival (and in the wave 
phasing) with respect to general relativity, similar to the dispersion of radio waves by interstellar plasma. 

Previous papers 0, 0, 0, IM El derived bounds on the graviton mass, on the Brans-Dicke parameter wbd and on 
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parameters describing more general scalar-tensor theories under the assumption that the compact objects do not 
carry spin. In this paper we investigate the effect of spin-orbit and spin-spin couplings both on the estimation of 
astrophysical parameters within general relativity, and on bounds that can be placed on alternative theories. We 
restrict our analysis to non-precessing spinning binaries, i.e. binaries whose spins are perpendicular to the orbital 
plane. The effect of non-aligned spins and the resulting precessions will be considered in future work. 

Within Einstein's general relativity, various authors have investigated the accuracy with which LISA can determine 
binary parameters including spin effects. Cutler 0] determined LISA's angular resolution and evaluated the errors 
on the binary masses and distance considering spins aligned or anti-aligned with the (orbital) angular momentum. 
Hughes investigated the accuracy with which the redshift can be estimated (if the cosmological parameters are 
derived independently), and considered the black-hole ring-down phase in addition to the inspiralling signal. Seto [s| 
included the effect of finite armlength (going beyond the long wavelength approximation) and found that the distance 
and angular resolution accuracy improve. This happens because the response of the instrument when the armlength 
is finite depends strongly on the location of the source, which is tightly correlated with the distance and the direction 
of the orbital angular momentum. Vecchio 0] provided the first estimate of parameters for precessing binaries when 
only one of the two supermassive black holes carries spin. He showed that modulational effects decorrelate the binary 
parameters to some extent, resulting in a better estimation of the parameters compared to the case when spins are 
aligned or antialigned with angular momentum. More recently, Hughes and Menou 11] studied a class of binaries, 
which they denoted golden binaries, for which the inspiral and ring-down phases could be observed with good enough 
precision to carry out valuable tests of strong-field gravity. 

These earlier works (except for JllJ) adopted analytical approximations to LISA's instrumental noise 0, augmented 
by an estimate of white-dwarf confusion noise |l2j in the low-frequency band. In this paper we model the LISA noise 
curve by a similar - albeit slightly updated - analytical approximation |13|. This noise curve has the advantage of 
being given in analytical form, and reproduces very well the salient features of numerical noise curves available online 
from the LISA Sensitivity Curve Generator (SCG) 0], a tool sponsored by the LISA International Science Team. 

Our central conclusions are as follows. Inclusion of non-precessing spin-orbit and spin-spin terms in the gravitational- 
wave phasing generally reduces the accuracy with which the parameters of the binary can be estimated. This is 
not surprising, since the parameters are highly correlated, and adding parameters effectively dilutes the available 
information. Such an effect has already been described within Einstein's general relativity in the context of ground- 
based detectors of the LIGO/VIRGO type [HHH- 

For example, for massive black-hole binaries at 3 Gpc, we find that 
including spin-orbit terms degrades the accuracy in measuring chirp mass by factors of order 10, and in measuring the 
reduced mass parameter by factors of order 20 - 100; while including spin-spin terms further degrades these accuracies 
by factors of order 3 and 5, respectively. For neutron stars inpiralling into intermediate-mass black holes (IMBH) 
with masses between 1000 and 10 4 solar masses, the corresponding reductions are factors of order 20 and 5 - 30 in 
chirp mass and reduced mass parameter, respectively, when spin-orbit is included, and additional factors of order 4 
and 7, respectively, when spin-spin terms are included. 

When we consider placing bounds on alternative theories of gravity, for technical reasons, we treat only spin-orbit 
terms. The source of choice to place bounds on the coupling parameter wbd of scalar-tensor gravity is the inspiral 
of a neutron star into an intermediate-mass black hole. We first reproduce results of earlier work apart from 
small differences arising from corrected normalization of the LISA noise curve. Including spin-orbit effects reduces 
the bound on wbd significantly, by factors of order 10 - 20. For example for a 1.4 Mq neutron star inspiralling into 
a 400 Mq black hole, the bound on cjbd goes from 8 x 10 5 to 40,000 when spin-orbit terms are included. The latter 
bound should be compared with the bound of 40,000 from Cassini measurements of the Shapiro time delay |l7|. 

The effect of including spin on bounding the graviton mass is more modest. In this case, the source of choice is the 
inspiral of binaries of massive black holes. For masses ranging from 10 5 M Q to 10 7 M Q , the reduction in the bound 
induced by the inclusion of spin-orbit terms is only a factor of 4 to 5. 

We consider the effect of spin terms on the angular and distance resolution of LISA. We find that spin couplings 
have a mild effect on the angular resolution, on the distance and, as a consequence, on the redshift determination 
for massive black-hole binaries. By contrast, for stellar mass objects inspiralling into intermediate-mass black holes, 
neither distance nor location on the sky is very well determined. 

LISA can observe massive black- hole binaries with large SNR out to large values of the cosmological redshift . If 
the corresponding mass and distance determinations are accurate enough, LISA will be an invaluable tool to study 
structure formation in the early Universe. Using Monte Carlo simulations we find that LISA can provide accurate 
distance determinations out to redshift z ~ 2 for source masses ~ 10 7 M Q , and out to z ~ 4 for source masses 
~ 10 6 M Q . Mass determinations strongly depend on an accurate treatment of spin effects. 

Finally, we study the effect of LISA's low- frequency sensitivity on the accuracy of estimating parameters for massive 
black-hole binaries (similar investigations can be found in [ljf and O). Below 10~ 4 Hz, the noise characteristics 
of LISA are uncertain. We show, however, that if LISA's noise spectral density can be uniformly extrapolated from 
10~ 4 to 10~ 5 Hz, then the accuracy of estimating both extrinsic parameters such as distance and sky position and 
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intrinsic parameters such as chirp mass and reduced mass, as well as the graviton mass, can be significantly higher, 
especially for higher-mass systems. 

The paper is organized as follows. In Sec. Ill Al we discuss the procedure for estimating binary parameters and the 
parameters of alternative theories when we average over all sky directions and binary orientations; this essentially 
ignores modulational effects due to the motion of the spacecraft. In Sec. lIIBl we relax this assumption, and present the 
relevant equations for estimation for a given source direction and orientation. In Sec. Ill CI we discuss the LISA noise 
curve to be used. Section IlIII presents our results. In Sec. IIII XI we show the results for estimates assuming averaging 
over directions. In Sec. IIII Bl we carry out a Monte Carlo analysis of 10 4 binaries distributed over the angles describing 
the relative orientation of the binary with respect to LISA and discuss the accuracy with which binary parameters 
can be estimated. Though somewhat more accurate, this procedure is still affected by various approximations. In 
Sec. IIII CI we use Monte Carlo simulations to investigate the dependence of parameter estimation on the redshift 
of the source. In Sec. IIII Dl we study the effect of the LISA low-frequency noise. Section llVl summarizes our main 
conclusions. In Appendix lAl we summarize for completeness the main equations used in this paper to describe LISA's 
configuration, orientation and response, as derived in Ref. @. In Appendix[5]we discuss some subtleties in estimating 
binary parameters within the Fisher matrix formalism in alternative theories of gravity when we include spin effects. 

Throughout this paper we use units in which G = c = 1 . 



II. ESTIMATION OF PARAMETERS IN NON-PRECESSING SPINNING COMPACT BINARIES 

We assume that two independent Michclson outputs h a (t) with a = I, II can be constructed from the readouts of 
the three LISA arms if the noise is totally symmetric (see, for example, pj). We take two approaches to estimating 
parameters. In the first approach (Sec. Ill Afl . we assume that the orientation and location of the source and the 
orientation of LISA are unimportant in estimating intrinsic parameters such as masses or theory-dependent parame- 
ters. These orientation dependences are contained in a number of angle-dependent functions, called pattern functions. 
Accordingly, we derive results using only one Michelson output, and we average over the pattern functions. 

In the second approach (Sec. Ill 5|) . we are also interested in the accuracy of determination of direction and distance 
to the source, and thus we do not wish to average a priori over pattern functions. Instead, we carry out Monte Carlo 
simulations of measurements using a population of sources across the sky, and we study the distribution of accuracies 
of parameter estimation. In this case, we use both one and two Michelson detectors. We use the by now standard 
machinery of parameter estimation in matched filtering for gravitational wave detection that has been developed by 
a number of authors [TEL I20I I2II.I2H . 



A. Parameter estimation using pattern-function averaging 



The Fourier transform of the waveform for one Michelson LISA detector, in the stationary phase approximation 
(SPA), and after averaging over the pattern functions, is given by 

M/) = y^/- 7/6 e W) , a -I, II, (2.1a) 
1 M 5/e 

A = = , — — , (2.1b) 

where / is the frequency of the gravitational waves, M = rf^M is the "chirp" mass, with M = mi + 1712 and 
i] = mi?ri2/Af 2 , and Z?l is the luminosity distance to the source. 

We have adopted the standard "restricted post-Newtonian approximation" for the waveform, in which the amplitude 
is expressed to the leading order in a post-Newtonian expansion (an expansion for slow-motion, weak-field systems 
in powers of v ~ (M/r) 1 / 2 ~ (irAif) 1 ^ 3 ), while the phasing ip(f), to which laser interferometers are most sensitive, 
is expressed to the highest post-Newtonian (PN) order reasonable for the problem at hand. For binaries with spins 
aligned (or anti-aligned) and normal to the orbital plane, this is a valid approximation because the amplitude varies 
slowly (on a radiation reaction timescale) compared to the orbital period. But when the spins are not aligned, 
modulations of the amplitude on a precession timescale must be included. Such modulations are beyond the scope of 
this paper. 

The phasing function ijj(f) is known for point masses up to 3.5 PN order |23l l24j . But spin terms are known only 
up to 2PN order, so to be reasonably consistent, we will include in the phasing point-mass terms only up to this same 
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2PN order. The needed expression for the phasing is 

+ (j§ + f<i) <,-*"(*Mff'°-m*v-"'(*Mf) + i0,,- 1 "(*Mf) 

The structure of the phasing function is as follows: the first two terms are related to the time t c and phase (f> c of 
coalescence; they are parameters that essentially establish where the waveform begins or ends. The prefactor of the 
expression in braces, together with the first term ("1") inside the braces, is the standard phasing from the lowest- 
order quadrupole approximation of general relativity. Inside the braces is a post-Newtonian expansion in powers of 
v <~ (irMf) 1 / 3 . The second term is the contribution of dipole gravitational radiation in Brans-Dicke theory. Let us 
define the scalar charge of the i— th body by a, = a&i — a(l — 2s%), where a 2 — 1/{2lobd + 3) ~ (2cjbd)~ 1 in the 
limit lobd 3> 1 ; and is called the sensitivity of the i— th body (a measure of the self-gravitational binding energy 
per unit mass). Then the coefficient in the dipole term is S = (a\ — &2)/2. The fact that it is dipole radiation 
means that it is proportional to v~ 2 compared to the quadrupole term, but the small size of S and the large current 
solar-system bound on wbd make this a small correction, nevertheless. The third term in the braces is the effect of a 
massive graviton, which alters the arrival time of waves of a given frequency, depending on the size of the graviton 
Compton wavelength X g and on a distance quantity D, defined below. The remaining terms in the braces are the 
standard general relativistic, post-Newtonian terms, including spin effects. 

The quantities and a represent spin-orbit and spin-spin contributions to the phasing, given by 



12 

i=l 



2 

Til 

113^ + 75, 



L • Si , (2.3a) 



a = ^XiX2 (-247Si • S 2 + 721L • SiL • S 2 ) , (2.3b) 

where Sj and L are unit vectors in the direction of the spins and of the orbital angular momentum, respectively, and 
Sj = XiTnjSi- For black holes, the dimcnsionless spin parameters \i must be smaller than unity, while for neutron 
stars, they are generally much smaller than unity. It follows that \0\ < 9.4 and \a\ < 2.5. 

We assume that any modifications to the post-Newtonian general relativistic terms listed in the phasing formula that 
might be generated in Brans-Dicke theory or in massive graviton theories will be of order 1/wbd <C 1 or l/(/A s ) <C 1 
relative to those terms, and hence we will ignore such corrections. 

In this paper we denote by M and M the observed chirp and total masses. They are related to masses measured 
in the source rest frame by 

M = (1 + z)M SOU rcc , M = (1 + z)M SOUIce , (2.4) 

where z is the cosmological rcdshift. 

Henceforth, to simplify the notation we define 

w = — , (2.5a) 

WBD 

To estimate the binary and gravitational theory parameters, we use the standard technique of parameter estimation in 
matched filtering. By maximizing the correlation between a template waveform that depends on a set of parameters 
9 a (for example, the chirp mass M.) and a measured signal, matched filtering provides a natural way to estimate the 
parameters of the signal and their errors. With a given noise spectral density for the detector, S n (f), one defines the 
inner product between two signals h\ (t) and h 2 (t) by 

(W^2jT ^jfV ( 2 - 6 ) 
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where h\(f) and h,2(f) are the Fourier transforms of the respective gravitational waveforms h(t). The signal-to-noise 
ratio (SNR) for a given h is given by 



p[h] = (h\h) 



1/2 



(2.7) 



If the waveforms may be characterized by a set of parameters 9 a , then one defines the "Fisher matrix" r o j, with 
components given by 



dh dh 
d6" 1 W 



(2.8) 



In the limit of large SNR, if the noise is stationary and Gaussian, the probability that the GW signal s(t) is charac- 
terized by a given set of values of the source parameters 9 a is 



p(6\s) =p(°)(0)exp 



-^T ab A9 a A9 b 



(2.9) 



where A9 a = 8 a — 9 a , and p(°\0) represents the distribution of prior information. An estimate of the rms error, A9 a , 
in measuring the parameter 9 a can then be calculated, in the limit of large SNR, by taking the square root of the 
diagonal elements of the inverse of the Fisher matrix, 

A6» a = Vs^, £ = r _1 . (2.10) 

The correlation coefficients between two parameters 6 a and 9 b are given by 

c ab = S afc /VS aa S fc6 . (2.11) 

We may wish to take into account our prior information on the maximum spin; we do this in a crude way by assuming 



p (o) (0) cx exp 



i(/3/9.4) 2 - i(a/2.5) 2 



The following derivatives of h will be needed: 

dh 



d\nA 
dh_ 

dt c 
dh 

d(j> c 
dh 
dw 
dh 

W 9 
dh 
d\nM 

dh 



d\nr/ 
dh 

dp 

dh 
cV 



h, 
2-nif h , 



-ih, 
5i 



3584 
i 

irMf 

5i 
128 



SW ,h {KMf)- 7 l*~h, 
h, 

{nMf)- 5/3 (K 4 v- 2 + 1 + A 4 v 2 + B 4 v 3 + C 4 v 4 ) h , 



-— (ttM/)- 5 / 3 (K 5 v- 2 + A 5 v 2 + B 5 v 3 + C 5 v 4 )h, 
9b 

~V- 3/5 (*Mf)- 2 /*h, 

-^n-^inMfy^h, 
04 



(2.12) 

(2.13a) 
(2.13b) 
(2.13c) 
(2.13d) 
(2.13e) 
(2.13f) 
(2.13g) 
(2.13h) 
(2.13i) 
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where here we denote v — (ttM/) 1 / 3 and 



A 4 

A 5 
B 5 



S 2 



w , 



12 
4 /743 
3 V 336 



11 

T 



128 



:(/3-47r), 

/ 3058673 5429 617 



V 1016064 
35 2 
56 

/743 
V 168 



1008 



n- 



144 



n 



■ zu , 



33 



27 

y(/?-47r), 



= 18 



/3058673 5429 617 



V 1016064 4032 



96 



rj — a 



(2.14a) 
(2.14b) 
(2.14c) 
(2.14d) 

(2.14e) 
(2.14f) 
(2.14g) 
(2.14h) 



For all integrals appearing in the Fisher matrix we will pick the final frequency, or the upper limit of integration to 
be /fi n = min(/isco ; /end)- Here /rsco is twice the conventional (Schwarzschild) frequency of the innermost stable 
circular orbit for a point mass, namely /isco = (6 3 / 2 7rM) -1 , and / cn d = 1 Hz is a conventional upper cutoff on the 
LISA noise curve. The initial frequency /j n in the integrals of the Fisher matrix is determined by assuming that we 
observe the inspiral over a time T bs before the ISCO, and by selecting a cutoff frequency below which the LISA noise 
curve is not well characterized. Our default cutoff is /i ow = 10~ 5 Hz; in Sec. IIII Dl we analyse the effects of increasing 
this cutoff frequency to reflect a less optimistic understanding of LISA's low frequency noise. The initial frequency is 
then given, in Hz, by the larger of these frequencies, 



fin = max /i ow , 4.149 x 10 5 



M 



l0 e Mr. 



-5/8 



lyr 



-3/8 



(2.15) 



The frequency at a given observation time is calculated using the quadrupole approximation for radiation damping. 
In our calculations we assume that T Q b s = 1 yr. 

Since we anticipate setting only lower bounds ohlobd and X g , we choose the nominal values zu — and [3 g = in Eqs. 
(|2.14|) . For simplicity, we will also assume that we are estimating spins in the case where spins are dynamically small. 
This is generally the case for neutron stars (see (2^ for discussion); for black holes, it means that we are considering 
only slowly rotating (non extremal) black holes. Consequently we also choose the nominal values f3 = a = in Eqs. 
(I2~T4T) . 



For a zero-spatial-curvature Universe (£l K — 0, Q\ 

1 + 



D, 



H 



Qm = l)j the luminosity distance is given by 
dz' 



[n M (i 



1/2 



The quantity D appearing in Eq. 12.5bf) is defined by (f2 K = 0) 

l + z f z dz' 



D 



o (1 + z') 2 [^ m (1 + ^') 3 + ^a] 



1/2 ' 



(2.16) 



(2.17) 



(see Eq. (2.5) of Q). For the Hubble constant we assume Hq = 72 km s _1 Mpc -1 , according to the present 
observational estimates |26| . 

A useful quantity to characterize the effect of the various terms (Brans-Dicke, massive graviton, spin couplings and 
PN corrections) on the evolution of the GW frequency is the number of GW cycles accumulated within a certain 
frequency band. This quantity is defined as: 



■A/gw = 



/in 



Idf. 
f 



(2.18) 
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To derive the number of cycles contributed by individual terms in the phasing, we use an expression for / that includes 
post-Newtonian GR terms, plus the Brans-Dicke and graviton-mass contributions, given by 



df 96 



dt 5ttM'- 



(irMf) 11/3 jl + r? 2 / 5 (nM fy 2 ^ + ^ (ttM/) 2 / 3 + PNcorrectionsj , (2.19) 



where the "PN corrections" up to 2PN order, including spin terms, can be found in [25|]. When we include the massive- 
graviton term, the frequency / and time t appearing in Eq. (|2.19l) should be considered as the arrival frequency and 
time, respectively (the number of gravitational-wave cycles due to the massive graviton being an effective number of 
cycles seen by the observer at the detector location). 

Contributions of individual terms in the integral l|2.18|) are generally considered significant if they exceed one wave 
cycle over the observation time. For various examples of the two source targets discussed in this paper, NS inspiral 
into IMBH, and MBH binaries, we show the individual contributions to the number of cycles, along with the initial 
and final frequencies, in Tables HI and ITT1 respectively. 

With the restricted post-Newtonian form for h in Eq. I|2.1|l . we can express the SNR J (p 2 ) in the form 



M\ 5/6 1/2 /lGpc\ //■/-* 3 /-*/3 



o-^j i"*(^)tL ikw*- (2 - 20) 

where angular braces mean that we are averaging over LISA pattern functions. 



B. Parameter estimation without averaging over pattern functions 

In this section we consider parameter estimation without averaging over the relative orientation of the binary with 
respect to LISA. We assume, as in [tJ, that two independent Michelson outputs can be constructed from the readouts 
of the three LISA arms if the noise is totally symmetric. The signal measured by LISA, h a (t) with a = I, II, can be 
written as: 

\/3 2m tyi f \ 

h a (t) = ^ cos (j o f^') dt ' + VvAt) + Vn(t)) > ( 2 - 21 ) 

where r(t) is the relative distance between the two compact bodies, (p p . a (t) is the waveform polarization phase [see 
Eq. ilA4a(l ] and fD(t) the Doppler phase [see Eq. (|A4b|l ]. A a (t) is defined by 



A a (t) = yf[l + (L • n)2]2 F + 2 + 4(L • n)2 F Q X 2 , (2.22) 

where L is the orbital angular momentum unit vector, and n is a unit vector in the direction of the source on the 
sky. The quantities F£' y are the pattern functions, defined by Eqs. I|A1(I and (|A2() . The Fourier transform of the 
measured signal can be evaluated in the stationary phase approximation, since A a (t), </? P)Q ,(t) and <po(t) vary on time 
scales on the order of 1 year (thus much larger than the binary orbital period ~ 2//) . The result is 

~ ha (f) = r v& e i * (/) ||i Q (t(/))| e - l (^w))+<p D w))) t (2 23) 

where to 2PN order (including also the Brans-Dicke parameter and the graviton-mass term) t(f) is given by 



S 2 w 



^ (.Ml)- 2 '* - {7rM f) 2 ^ 



+2 



/ 3058673 5429 617 2 \ _ 4/5 . iJX , 4/3 
+ 1008 ^+144^ ' ( ^ /) 7 



V 1016064 



(2.24) 



In Appendix [S] using equations of Ref. 0, we show how to express the angular parts of h a (f) in terms of the angles 
9si 4>Si ®l, 4>l, which describe the source location and orbital angular momentum direction in the reference frame 
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attached to the solar system barycenter. To evaluate the Fisher matrix we use the derivatives with respect to the 
parameters Ai, n, 0, a, <p c , t c , hi A, w and g , given by Eqs. (|2.13|) . We also determine analytically the angular 
derivatives with respect to 8$, <fis, 0l, 4>l using formulas given in Appendix lAl The final results are lengthy and 
unenlightcning, so we do not write them down here. We choose to evaluate the angular derivatives analytically since 
this is likely to be more accurate than the numerical finite-differencing adopted in Refs. [7j,|8|]. As before, we choose 
the nominal values vo — (3 g — (3 = a = 0. 
The non-averaged SNR is (a = I, II) 



p a {0s,4>s,h,4>L) = 7.807 xHT 23 (iQ 5/ V/2 {^>j ^ /C " d lAl{t{J);h, 4>s, h, ^)^~df ■ (2.25) 

For some binary-mass configurations we estimate the parameters using the two LISA detectors. In this case the 
Fisher matrix is 

r tot _ ( dhi dhA ( dh u dh u \ 
and the rms error, A9 a , in measuring the parameter 8 a is A9 a = VS aa with S = [r tot ] 1 . The total SNR is 

Ptat = \Jp\ + Pu- 

We expect that estimates obtained for parameters such as M, n, (3, cr, w and f3 g when we do not pattern average 
will not differ qualitatively from those obtained using pattern averaged templates. As we will see, those parameters, 
which appear in the phasing of the signal, are relatively uncorrelated with the parameters appearing in the amplitude, 
such as 6s, <j>s, $l- 



C. Noise curve for the LISA instrument 



The non-sky-averaged noise spectral density of LISA depends on the relative orientation between the instrument and 
the source, and it is very hard to implement in estimating binary parameters. Generally, the LISA community has 
been using the so-called sky-averaged spectral density S^ A [see e.g., Ref. [2jJ and the LISA Pre-Phase A Report]. 
The sky-averaged spectral density is computed by a combination of three factors, including: (i) the raw spectral noise 
density S n , (ii) the gravitational- wave transfer (response) function R and (iii) the noise transfer (response) function 
R n . They combine together in [2^ 

S SA = ^Rn (22?) 

In this paper we are also interested in determining binary parameters without averaging over the source location, so 
we are not allowed, in principle, to use 5*^ A . To overcome this difficulty we evaluate an effective non-sky-averaged 
spectral density which gives the correct result at low frequency, but is only approximately valid in the high-frequency 
region. In the low frequency limit, the GW transfer function used in the LISA Sensitivity Curve Generator [lj] is 
R = 4(V3/2) 2 l/5 = 3/5, where the factor (V3/2) 2 comes from the LISA arms being at 60°, the factor 1/5 is due to 
the sky-average of the pattern functions ((Ff x ) = 1/5) and the factor 4 depends on the particular read-out variable 
used. Since our definition of the GW signal already includes the factor y/3/2 [see Eqs. Ij2.1a|l and l)2.23[l ]. to obtain 
the effective non-sky-averaged spectral density we must multiply A by (\/3/2) 2 /5 = 3/20. The final result is: 



f \ / f 



9.18 x 10~ DZ -4— + 1.59 x 10~ 41 + 9.18 x 10 



1 Hz J V 1 Hz 



Hz" 1 , (2.28) 



and has been obtained also in Ref. [lj|. We estimate white-dwarf confusion noise following [13j . which uses results 
from |29l l30j | : the galactic contribution is approximated as 



,sr l m=2.1xlO- 45 ( T ^) Hz" 1 , (2.2!)) 



and the contribution from extra-galactic white dwarfs as 



-7/3 



S e *- Sa \f) = 4-2 X 10- 47 ( JL- ) Hz' 1 . (2.30) 
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We compute the total (instrumental plus confusion) noise as 

S h (f) = min {5, NSA (/)/exp (- nTj ssion dN / df) , S^ A (f) + Sf(f)} + S^if) • (2-31) 

Here dN/df is the number density of galactic white-dwarf binaries per unit gravitational-wave frequency, for which 
we adopt the estimate 

dN „ „ _ n (\ Hz^ 11/3 



i/f 2 x 10- 3 Hz" 1 ; (2.32) 

A/ = T~ ission is the bin size of the discretely Fourier transformed data for a LISA mission lasting a time T m i ss i on ; 
and k ~ 4.5 is the average number of frequency bins that are lost when each galactic binary is fitted out. The factor 
exp ( — KT~^„- n „dN / df) thus represents the fraction of "uncorrupted" bins where instrumental noise still dominates. 
At variance with |l3j . in our calculations we always assume that the duration of the LISA mission T m i ss i on = 1 yr, 
consistently with the choice we made for the observation time T b s in Eq. (|2.15l) . The analytic root noise spectral 
density curve (|2.31|) used in this paper is shown in Fig . HI together with the corresponding root noise spectral density 
curve from the LISA Sensitivity Curve Generator [lj. The SCG curve shown is obtained using the nominal values 
SNR=1, arm length = 5 x 10 9 m, telescope diameter = 0.3 m, laser wavelength = 1064 nanometers, laser power 
= 1.0 Watts, optical train efficiency = 0.3, acceleration noise = 3 x 10~ 15 ms~ 2 Hz -1 / 2 , and position noise budget 
= 2 x 10~ n mHz" 1 ' 2 , with position noise setting the floor at high frequency. The data returned by the SCG is then 
multiplied by a/3/20 to obtain the effective non-sky averaged curve shown in Fig. ^ 



III. RESULTS OF PARAMETER ESTIMATION 

A. Estimates using pattern-averaged templates 

We begin with neutron-star inspirals into intermediate mass black holes. These are the best sources for bounding 
scalar-tensor gravity, for the following reasons. In scalar-tensor theory, dipole gravitational radiation is controlled by 
the difference S = (&i — &2)/2 in the rescaled scalar charge d, between the two bodies. We recall that on = (1 — 2sj) 
and 

_ / <9(ln TTjj ) \ 

Si -{d(lnG eB )) N ' {3A} 

where mi is the total mass of the body, G e ff is the effective gravitational constant at the location of the body (which 
is related to the value there of the Brans-Dicke scalar field) and the subscript N denotes that the number of baryons 
is held fixed. For neutron stars, Sj can be substantial, and thus on can differ markedly from unity (pn ~ 0.6 — 0.8), 
but because it is only weakly dependent on the NS equation of state and mass |3l| , the difference S for NS binaries 
is typically S < 0.05, so NS-NS binaries do not provide interesting bounds on Brans-Dicke theory Q (see however 
for discussion of more general scalar-tensor theories). Because of the no-hair theorem, for black holes «bh = 0, so 
BH-BH binaries cannot be used to put bounds on the Brans-Dicke parameter via dipole radiation. Therefore, following 
previous papers on the subject we only consider NS-BH binaries as sources for this purpose. Furthermore, 

as shown earlier 0, Q, inspiral into lower-mass black holes gives the most promising bounds, primarily because 
more cycles are observed in a given period of integration in that case. The event rate of such inspirals involving 
intermediate-mass black holes is uncertain, but is likely to be very small [32L l33|: only a lucky detection of such an 
inspiral will lead to a suitable test. White dwarf-BH binaries could also be used to test Brans-Dicke theory, since 
Oiwo ~ 1 {swd ~ 10~ 4 ), so that 5 ~ 0.5, except that tidal effects will play a role in the late stages of the inspiral, 
depending on the mass of the black hole (for discussion, see Q)- 

For concreteness we focus on four NS-BH binaries, setting the NS mass Mns = 1.4M© and considering black holes 
of mass M BH = 400, 1000, 5000 and 10 4 M Q . The borderline between massive and supermassive BHs is hazy, but 
we choose not to consider NSs inspiralling into "supermassive" BHs with M = 10 5 -10 8 Mq. Even in the context of 
pure general relativity, our approximation that the binary orbits be circular is expected to be unreliable for these 
high-mass cases: for high mass ratios the binary is likely to be formed by capture of the smaller body into the larger 
one, and the eccentricity will not be washed out by radiation reaction. Adding eccentricity complicates the analysis 
to a level that is beyond the scope of this paper, and we plan to return to this problem in the future. 

We first consider the inspiral of these four representative binaries within general relativity (omitting the BD term 
in the phase). From the initial and ending frequencies listed in Tabled we see that these binaries sweep through 
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the high frequency part of the LISA band, say from ~ 10~ 2 Hz up to ~ 1 Hz. In Table [TTT1 we list the errors and 
correlation coefhcients that are obtained when we truncate the phasing formula at various PN orders and include 
spin-orbit and spin-spin effects. For consistency, at 1PN order we do not include spin effects, at 1.5PN order we 
include only spin-orbit effects, and at 2PN order we include both spin-orbit and spin-spin effects. All results assume 
p = 10. In the left panel of Fig. [21 we show the corresponding luminosity distance as a function of the black hole mass 
M BH . 

From Table lTTTl we see that the errors on all parameters increase considerably when spin effects are taken into account. 
This applies in particular to the chirp mass M and the parameter r/. This spin-induced degradation in parameter 
estimation has long been known it occurs because (in the absence of precessional effects) the parameters are 
highly correlated, so that adding parameters effectively dilutes the available information. 

For technical reasons, when we consider alternative theories of gravity we only include spin-orbit effects. If in 
addition we include spin-spin effects the dimensionality of the Fisher matrix increases, and the matrix inversion 
required to obtain the correlation matrix appears to be unreliable. This issue is addressed in Appendix iBl 

In Table ITVl wc show errors and correlation coefficients for NS-BH binaries at 2PN order when we include the BD 
term. For nonspinning binaries the results are similar to Table I in Ref. pj , except that those authors used templates 
at 1.5PN order and did not take into account the factor v3/2 which appears in Eq. p.lajl . The BD term is highly 
correlated with M. and r\ (c M ™ and d 1 ^ are both quite large). Correspondingly, the error on both M. and ij increases 
by roughly one order of magnitude with respect to the "general relativistic" values listed in Table ITTT1 We also compute 
the BD bound obtained by inverting only the diagonal element of the Fisher matrix. This "uncorrelated" bound 
WBD,unc is always about two orders of magnitude larger than the actual value we obtain by inverting the full Fisher 
matrix in the absence of spins. 

Notice also that the BD bound decreases with increasing black hole mass. This can be partially understood by the 
following argument: the derivative of the GW signal with respect to w is proportional to A/^g Mg^ 3 [see Eq. (|2.13d|) 
and use Mbh 3> Mns]- Therefore the derivative decreases as Mbh increases; the integration over the frequency range 
(which also depends on mass) modifies this dependence somewhat, but the final conclusion is that the higher the BH 
mass, the lower the Brans-Dicke bound. 

From Table IIVI we also see that non-precessional spin effects reduce considerably the bound on the Brans-Dicke 
parameter. For example, for a (1.4+1000) Mq binary the bound decreases by a factor 10 (from ~ 2x 10 5 to ~ 2 x 10 4 ) 
when we include the spin-orbit term. We will see later ( Table IVTT)l that a further reduction of a factor ~ 2 comes 
from inclusion of effects related to the orbital motion of LISA. We also found that including priors, that is, assuming 
that we know a priori from general relativity that the compact objects' spins are bounded from above (see [l^ for a 
discussion) has completely negligible effects when we include spin-orbit terms. 

In setting bounds on the graviton mass, we consider massive and supermassive binary black holes with M = 
10 4 -10 7 A/©. The derivative of the GW signal with respect to (3 g is inversely proportional to the chirp mass A4, so for 
comparable-mass binaries, the higher the total mass the higher the graviton-mass bound. However, as seen from the 
initial and ending frequencies in Table ITT1 for Mbh > 10 7 Mq, the binary sweeps through the low- frequency end of 
the LISA band below 10 -4 Hz, where the predicted sensitivity of LISA is not very robust at present. The sensitivity 
in this low-frequency regime depends on how efficiently the acceleration noise can be reduced. In Sec. Ill Bl we shall 
investigate the effect on the estimation of the parameters and on the graviton-mass bound, if the LISA noise curve 
can be trusted only down to a lower frequency /i ow ~ 10~ 4 or 5 x 10 -5 Hz. 

In Tables and IVII we list the errors and correlation coefficients when binaries with high, comparable masses 
are detected using pattern-averaged templates at 2PN order. Table |V] shows results for pure general relativity, with 
spin-orbit and spin-spin effects included. Table IVll shows results when a massive graviton term and a spin-orbit term 
are included. As in the BD case, we do not show results for a massive graviton combined with spin-orbit and spin-spin 
effects, because the inversion of the large Fisher matrix in this case appears to be unreliable (see Appendix [HJ . 

As in the case of NS-BH binaries, adding new parameters causes a degradation in the accuracy with which we 
can estimate parameters. All the values reported in Tables [3 and IVII are obtained for binaries at 3 Gpc. The 
corresponding SNR for equal mass BH-BH systems is shown in the right panel of Fig. [2] as a function of the total 
mass of the binary. For M > 1O 6 M0 we observe the appearance of a relative minimum, corresponding to the range 
of frequencies in which white-dwarf confusion noise dominates over instrumental noise. 

Although we only report results for the currently favoured values of the cosmological parameters, we verified 
that the upper bound on the graviton wavelength depends only weakly on the underlying cosmological model. The 
"uncorrelated" bound A 9iUnc obtained by inverting only the diagonal element Tp g p g of the Fisher matrix is about 
one order of magnitude larger than the result obtained by inverting the full Fisher matrix (for wbd the difference 
was about two orders of magnitude). Compared with the case of scalar-tensor theories, bounds on massive graviton 
theories seem to be less sensitive to correlations among different parameters. 
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B. Estimates using templates without pattern averaging 

To assess the effect of pattern-averaging on parameter estimation, and also to determine how accurately LISA can 
measure source locations and luminosity distances, we adopt the non-averaged templates of Sec. Ill Bl and perform 
Monte Carlo simulations using a population of sources across the sky. We consider in detail two systems: (i) a NS-BH 
binary with mass of (1.4 + 1000)Af Q observed with a single-detector SNR pj = 10 (a typical target system used to 
place bounds on the BD parameter), and (ii) a BH-BH binary with mass of (10 6 + 10 6 )M Q at distance Dl = 3 Gpc 
(a typical system in the context of massive graviton theories). 

For each of these systems we distribute 10 4 sources over sky position and orientation. We randomly generate the 
angles <f>s, 4>l in the range [0, 2tt] and ps = cos#s, Pl — cos9l in the range [—1,1]. As in to generate random 
numbers we use the routine RAN2 |34j |. Computing and inverting the Fisher matrix for 10 4 binaries typically takes 
^8—15 minutes (depending on the dimensionality of the matrix) on an ordinary laptop. This is much faster (by 
a factor ~ 500) than previous Monte Carlo simulations [3{|. A marginal difference with previous codes is that we 
compute angular derivatives analytically instead of using finite differencing, but the major improvement is due to 
our use of a numerical integrator based on spectral methods, namely the Gauss-Legendre routine GAULEG [34| . 
Numerical experiments show that ~ 600 points in the spectral expansion are sufficient to obtain an accuracy of one 
part in 10 4 in all parameter errors. This is true even when we use splines to interpolate tabulated data from the 
LISA Sensitivity Curve Generator, instead of adopting the analytical noise curve of Eq. I|2.31[) . When the waveform 
contains a large number of highly correlated parameters, computing the inverse of the Fisher matrix can be numerically 
difficult. The method we used to check the robustness of our results is described in Appendix 151 

Once we have computed the errors for all 10 4 binaries we group them into Ab; ns bins depending on the (logarithmic) 
distribution of their errors: a binary belongs to the j-th bin if the error on some parameter X satisfies 

I JVbins J I JVbins J 

for j = 1, . . . , Abi ns . In this paper we fix iVbms = O nce we have binned the data, we normalize the binaries in 
each bin to the total number of binaries to get a "probability distribution" of the error on the variable X . 

In Figs. |31 and 01 we show the resulting histograms for a NS-BH binary of (1.4 + 10 3 )M Q with pj — 10. The 
plots contain various histograms for parameter estimations made when spins are absent, when spin-orbit is included, 
and when both spin-orbit and spin-spin are included. The histograms come in pairs: in each case the solid-line 
histogram refers to measurements carried out with only one data stream from the Michelson interferometer I, the 
dashed histogram refers to measurements made when both data streams from Michelson interferometers I and II are 
combined. Not surprisingly, the accuracy is improved with the use of two outputs, very roughly by a factor of order 
V2 in most cases. 

In Fig. |2| we plot the probability distribution for the luminosity distance ADl/Dl and for the angular resolution 
AOs, defined as 

AO s - 2n\ Binflsl {z Sg g g V$ g } a ~ ^s}^ ■ ( 3 - 3 ) 

In Fig. 0| we plot the distributions for the chirp mass AAA/A4, the reduced mass Ap/p, the spin parameters A/? 
and Act, and the bound on the BD parameter cjbd- The distribution for Ap/p can be obtained from the errors on r\ 
and M. by error propagation, taking into account that the correlation between the two mass parameters, as defined 
in Eq. (|2.11() . can be large: 

Ap 

The plots contain histograms evaluated when spins are absent, when spin-orbit alone is considered and when both 
spin-orbit and spin-spin are included. We complement these plots by Table IVIII which shows the average errors 
obtained summing over all binaries, both with and without the BD term. Notice that this procedure is different 
from averaging over the sky without taking into account the orbital motion of LISA. For each model, the first line 
in Table IVTT1 refers to errors obtained averaging over all binaries and using only detector I, the second line refers to 
an average over all binaries using both detectors, while the third line reproduces, for comparison, the corresponding 
pattern-averaged results from Tables llTfl and llVl In general, the pattern-averaged procedure gives good qualitative, 
albeit systematically low estimates of the measurement errors, compared to the Monte Carlo results. For stellar mass 
inspirals into intermediate-mass black holes, both the angular resolution and the distance determination accuracy are 
poor. The Monte Carlo simulation gives rather broad probability distributions, shown in Fig. El with minimum error 
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in distance around 0.1, but with a tail extending up to ADl/Dl ~ 10 or so. The values for Ails in steradians look 
rather small, but when expressed in arcminutes, with 60s ~ [Ails x (3283/str)] 1 / 2 x 60 arcmin, they are substantial. 
For comparison, the angular diameter of the Moon (and the Sun) as seen from the Earth is ~ 30 arcmins. The angular 
resolution is degraded when spin terms are included, as is apparent in the top panels of Fig. El while the distance 
determination is relatively insensitive to the inclusion of spin terms. As first noticed by Cutler |7j, both the angular 
distribution and the distance determination improve when we use both detectors (dashed lines) instead of a single 
detector (solid lines). 

A noteworthy feature of the histograms in Fig. 0] is that the errors on M,, fi, /3, cr and wbd show a peculiar "double- 
peak" structure which is absent for the high-mass BH-BH binaries. We see this double-peak structure for the first 
time because our fast spectral integrator allows us to simulate a sufficiently high number of binaries, but we have 
no analytical understanding of this behavior. The inclusion of each spin-coupling term degrades the determination 
of both M. and ji by roughly one order of magnitude. The large reduction on the Brans-Dicke bound caused by the 
inclusion of non-precessing spins is one of the main results of this work (bottom panel of Fig. EJ- 

In Figs. Eland El we show similar histograms for a binary of total observed mass (10 6 + 10 6 )M© at a fixed distance 
of 3 Gpc, and in Table IVIIII we display average errors obtained by summing over all binaries. Supermassive BH-BH 
binaries can be observed at much higher redshifts than NS-BH binaries. LISA's accuracy in measuring the luminosity 
distance Dl can thus be exploited to infer the redshift z of the source, disentangling the mass-redshift degeneracy 
of the waveforms and allowing the determination of the masses in the source rest frame. Indeed, Hughes proposed 
to use gravitational wave observations in this way to map the merger history of supermassive black holes 8]. Once 
we have ADl/Dl we can compute the error on the binary's redshift by the following procedure. If we assume that 
the cosmological parameters il\ and H are known with an accuracy ~ 10%, we can use error propagation to get Az 
from ADl following Hughes @|: 



Az 
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dD L 
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(3.5) 



From Eq. (|2.16|) . and assuming that VLm + Ha = 1, we find the derivatives 
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(3.7) 



Thus Az is completely determined once we fix z, ADl/Dl, AH /H , Ha, AQa- We also compute the best possible 
redshift determination (Az)bcst that LISA could achieve assuming (perhaps not too optimistically) that by the time 
LISA flies the cosmological parameters are known to much better precision than LISA's distance determinations 
(Afl A ps , AH ps in Eq. 

Figure [SI shows the resulting probability distributions for the SNR, the luminosity distance ADl/Dl, the redshift 
determinations Az/z and (Az/z)bcst, and the angular resolution AQ$ m steradians. Unlike the NS-BH systems 
considered in Fig. the probability distribution of all these quantities (at fixed distance) depends very weakly on 
whether we include or omit spins, so we display only the spinless results. The SNR distribution has the same shape 
and average value for the two Michelson detectors; it increases (on average) by a factor ~ y2 when we use both 
detectors. Not surprisingly, the distribution of (Az/z)best and ADl/Dl are identical, apart from the scale factor 
from dD^/dz in Eq. I|3.6a|l . The shape of the distribution of ADl/Dl is different from that shown in Fig. El The 
distribution of Az/z is dominated by the 10 per cent errors assumed for the cosmological parameters, and shows only 
small effects of the distribution of ADl/Dl. At this relatively small redshift, the present uncertainty on cosmological 
parameters dominates over the accuracy of LISA 8]. For larger values of the redshift LISA distance and redshift 
measurements become less accurate, eventually dominating the error at some critical value of z that depends on the 
MBH masses (see Sec. lHIC)) . 
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Figure [5] shows distributions for the chirp mass AA4/A4, the reduced mass A/x//i, the spin parameters, A/3 and 
Act, in general relativity, and finally the bound on the graviton Compton wavelength when that term is included. 
The corresponding errors are listed in Table IVlTll It turns out in this case that the "pattern-averaged" approach of 
Sec. Ill Al provides very good estimates of the errors on the parameters M, rj, \ g , (3 and a (in fact, these estimates 
are almost identical to those obtained from the Monte Carlo simulations using both detectors, as can be seen by 
comparing the relevant values in Tables Ivl and IVll with Table IVriT|) . Quite remarkably, the errors on Dl, -M, f2s and 
z in the case GR+SO+SS are exactly the same as the corresponding errors for the case MG+SO, while the errors on 
r\ and /3 differ (a similar consideration applies to the angle- average calculation of Sec. IIII A"|l . 

The determination of M. and fj, is affected by the use of one or two detectors and by the inclusion of the spin 
couplings in the same way. Notice however that A/i//i is always about two orders of magnitude larger than AM/Ai. 
The determination of the spin-spin parameter a without large errors is only made possible by the use of two detectors. 

When wc consider alternative theories of gravity (either with a BD or a massive graviton term) we do not include 
spin-spin effects because the Fisher matrix becomes ill-conditioned and non-invertible (see Appendix [BJ) . The non- 
invertibility of the Fisher matrix in the case MG+SO+SS is probably related to the "degeneracy" of the errors in the 
cases GR+SO+SS and MG+SO. 

C. Redshift dependence of parameter estimation 

Inspiralling MBH can be observed by LISA out to enormous distances. If their masses and luminosity distances 
are determined with sufficient accuracy, LISA can be a source of information on the growth of structures at high 
redshift. In this context the redshift z can be large, and the distinction between observed masses and masses as 
measured in the source rest frame - given by the simple rescaling l|2.4[l - is important. LISA can only measure 
redshifted combinations of the intrinsic source parameters (masses and spins), so it cannot measure the redshift z. 
If cosmological parameters are known, Eq. (|2.16|) can be inverted to yield z as a function of Dl and LISA 
measurements of the luminosity distance can be used to obtain black hole masses as a function of redshift, thus 
constraining hierarchical merger scenarios. Alternatively: if we can obtain the binary's redshift by some other means, 
e.g., from an electromagnetic counterpart, then LISA measurements of Dl can be used to improve our knowledge of 
the cosmological parameters [3(| 153, H|| . 

These exciting applications depend, of course, on LISA's measurement accuracy at large redshifts. In the following 
we look at the redshift dependence of measurement errors for two representative MBH binaries having masses (10 6 + 
10 6 )-Mq and (10 7 + 10 7 )Mq as measured in the source rest frame. (This choice is at variance with the rest of the 
paper, where we fix instead the values of the measured masses at the detector.) We consider a zero-curvature-universe 
with (Ojvf = 0.3, fl\ = 0.7) and Ho = 72 km s _1 Mpc -1 , according to the present observational estimates. We also 
assume that the LISA noise can be extrapolated down to /i ow = 10~ 5 Hz; more conservative assumptions on /i ow could 
significantly affect our conclusions (see [lSl Il9l| and Sec. IIII 01 . We compute the errors, as in Sec. IIII 51 performing 
Monte Carlo simulations of 10 4 binaries for different values of the redshift and then averaging over all binaries. 

Fig. [7]shows the redshift dependence of the average errors on various quantities (all errors are computed using two 
detectors). The left panel corresponds to the general relativistic inspiral of a nonspinning binary. In the right panel 
we include, in addition, SO and SS terms. Solid (dashed) lines refer to MBH binaries having mass (10 6 + 10 6 )Af Q 
[(10 7 + 10 7 )M Q , respectively] as measured in the source rest frame. 

As expected from the discussion in Sec. IIII 51 distance determination and angular resolution are essentially indepen- 
dent of the inclusion of spin terms. The relative error on Dl for the lower-mass binary system is ~ 2% at z = 1, ~ 5% 
at z — 2 and ~ 11% at z = 4. This reduction in accuracy is due to the fact that the signal spends less and less time 
in band as the redshift is increased. We only consider values of the redshift such that the binary spends at least one 
month in band before coalescing: for the case (10 7 + 10 7 )M Q , this corresponds to z ~ 4. The distance determination 
error for this high-mass binary grows quite rapidly, being ~ 2% at z = 1, ~ 6% at z = 2 and ~ 21% at z = 4. LISA's 
angular resolution is rather poor even at small redshifts, and it rapidly degrades for sources located farther away, 
the degradation being more pronounced for higher-mass binaries. Better distance determinations can be obtained if 
we are lucky enough to locate the source in the sky by some other means: for example, associating the gravitational 
wave event with an electromagnetic counterpart. In this case ang les and distance would be decorrelated, allowing 
order of magnitude improvements in the determination of Dl |18j . We should also recall that in our discussion we 
are quoting average errors. Since the logarithmic distribution of ADl/Dl has a rather long tail at large values of the 
error (Fig. [SJ, distance errors in a specific detection could actually be much smaller than the average. 

Unimportant as they are for distance determination and angular resolution, spin effects have a dramatic impact 
on mass measurement accuracy. For our low-mass system, in the absence of spins the chirp mass can be measured 
with fantastic accuracy up to z = 10, the largest error being ~ 0.06%. Even including SO and SS effects and ignoring 
precession, the error on M. is only ~ 2.5% at z = 10. Errors are predictably larger for the (10 7 + 10 7 )M Q binary. 
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When we omit spin effects Ad can still be measured with an accuracy better than a percent out to z = 4, but when 
we include SO and SS terms the error at z — 2 is already ~ 6%. Our ability to measure the mass of both black holes 
is severely limited by the error on the reduced mass fi, which is always about two orders of magnitude larger than the 
error on Ad. Errors on fi for nonspinning binaries of (10 6 + 10 6 ) Mq are remarkably small if we ignore spin effects: at 
z = 10 the error is only 0.6%. Including both SO and SS terms things get much worse, and even at z = 1 the reduced 
mass error is ~ 6%. 

For (10 6 + 1O 6 )M0 general relativistic nonspinning binaries, a least-square fit of mass and distance errors in the 
interval z G [1, 10] yields: 

AM/M = (-1.1476 + 7.2356z + 5. 7376z 2 ) x 10" 6 , 

A/i/n = (-0.61431 + 1.9018z + 0.43721z 2 ) x 10" 4 , 
AD L /D L = (-0.65651 + 2. 6935z + 0.061595z 2 ) x 10~ 2 . 

It is important to remark here that in our study we are ignoring precessional effects. These effects induce modula- 
tions in the waveform, possibly improving the mass measurements in a significant way [rJ. The study of precession 
is therefore crucial to assess LISA's ability to measure MBH masses in galactic mergers. Such a study is beyond the 
scope of this paper. 



D. Effect of LISA's low frequency sensitivity on parameter estimation 

High-mass binaries sweep through the low-frequency region of the LISA band, where the LISA sensitivity will 
ultimately depend on design choices for the acceleration noise. To explore the possible consequences of such design 
choices, we studied the effects on the accuracy of parameter estimation if LISA were completely blind below some 
cutoff frequency f\ ow . I n ah our analyses to this point, we chose the default value for f\ ow to be 10 -5 Hz, but we now 
consider higher cutoff frequencies 5 x 10 -5 Hz and 10 -4 Hz. 

Table IVnTl shows that if /i ow = 5 x 10~ 5 Hz , parameter estimation for binaries of (10 6 + 10 6 ) Mq is essentially 
unaffected. But if /i ow = 10~ 4 Hz, the accuracies of estimating all quantities are degraded by factors between two 
and six, both for GR and for massive graviton theories. 

Since binaries of larger mass sweep through a lower frequency band, we expect the degradation to be even worse for 
such binaries. To be more quantitative, we analysed what happens when we increase the LISA cutoff frequency /i ow 
from 10 -5 Hz to 10~ 4 Hz for massive black hole binaries having total mass larger than ~ 10 5 Mq. Similar preliminary 
investigations of the effect of low-frequency LISA noise on the distance determination of massive black holes can be 
found in [l^. I n this subsection we consider for concreteness nonspinning, equal-mass binaries. All quoted values for 
the average errors have been obtained using Monte Carlo simulations with two detectors. 

The results are summarized in Fig. [SJ where we show errors computed both in general relativity (left panel) and 
in massive graviton theories (right panel). In the plots we display results for cutoff frequencies /i ow = 10 -5 Hz and 
/ low = 10- 4 Hz. 

Consider first the distance determination accuracy ADl/Dl one could achieve. For general relativistic binaries, 
as we move the low-frequency cutoff /i ow , ADl/Dl grows from 1.1 % to 1.8 % for a binary of (10 6 + 10 6 ) Mq. If 
we double the mass of each black hole, the corresponding increase is roughly twice as large - from 1.0 % to 2.6 %. 
For larger values of the mass and /i ow = 10 -4 Hz the binary does not spend much time in band, and the inversion 
of the Fisher matrix becomes problematic (see Appendix [BJl . The angular resolution is even more sensitive to the 
low-frequency cutoff. As the cutoff /i ow goes from 10 -5 Hz to 10 -4 Hz, Afig (in steradians) goes from 9.7 x 10~ 5 
to 4.8 x 10~ 4 for a binary of (10 6 + 10 6 ) Mq. If we double the mass of each black hole ACls correspondingly goes 
from 8.4 x 10~ 5 to 1.8 x 10~ 3 , increasing by a factor ~ 10 2 . We note another interesting feature emerging from 
Fig. 00 If present design choices allow LISA to be sensitive down to /i ow = 10~ 5 Hz, the distance error and angular 
resolution (at least for Dl = 3 Gpc) will be decreasing functions of M in the supermassive black hole mass range 
M € [10 6 — 10 7 ] Mq. This feature could be used to study the merger history of black holes at galactic centers, and 
to map structure formation in the early universe. 

Mass determinations are also strongly affected by low- frequency sensitivity, as observed in Refs. 0,0- For a 
(10 6 + 10 6 ) Mq binary in general relativity, a lower cutoff at /i ow = 10~ 4 Hz increases the error on the chirp mass 
by a factor 3 for a (10 6 + 10 6 ) Mq binary, while the error in /i correspondingly increases from 0.0066 % to 0.012 %. 
An error of the order of a percent on mass determination during coalescence could be a problem for the identification 
of "golden binaries" (binaries for which we can measure the total mass-energy lost to GWs 11]) with LISA. If in 
addition to a cutoff at /i ow = 10~ 4 Hz we also include a massive graviton term, the errors for binaries of total mass 
10 7 Mq become unacceptably large: ~ 0.6 % for the chirp mass and ~ 4 % for the reduced mass. 
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Next we consider how the low-frequency cutoff affects bounds on the graviton mass. The results are shown in Fig. El 
It turns out that the bound on the graviton Compton wavelength is affected by low frequency noise in the same way 
as the accuracy in distance and redshift determination . As we increase /i ow from 10 -5 Hz to 10~ 4 Hz, the bound on 
X g (in units of 10 15 km) drops from 49.4 to 29.5 for a binary of (10 6 + 10 6 ) Mq. For an equal-mass BH-BH binary of 
total mass 10 7 M© the corresponding reduction is from 67.9 to 10.3. Notice also that when we pick /i ow = 10 -4 Hz 
the bound on A g has a maximum for equal-mass binaries of total mass ~ 10 6 Mq. Being sensitive below 10~ 4 Hz is 
therefore important to put bounds on the graviton mass through observations of binaries more massive than this. 

In conclusion, design choices for the LISA low-frequency acceleration noise will have a dramatic impact on our ability 
to i) locate massive black hole binaries in the sky, ii) measure their masses, iii) use them as standard cosmological 
candles, and iv) bound the mass of the graviton. 

IV. CONCLUSIONS 

In this paper we analysed how the inclusion of spin couplings affects parameter estimation in observations of binary 
coalescence, and the bounds that can be placed on alternative theories of gravity, such as the Brans-Dicke and massive 
graviton theories. Extending previous investigations 0, EH 0, 01 > we also took into account the dependence on the 
four angles describing the source location and the direction of the (orbital) angular momentum, performing large-scale 
Monte Carlo simulation of 10 4 binaries. 

We found that the bound on the Brans-Dicke parameter (and therefore also the bound on parameters describing 
more general scalar-tensor theories, such as those considered in Ref. 0) is significantly reduced by spin-orbit and 
spin-spin couplings, while the bound on the graviton Compton wavelength is only marginally reduced. As expected, 
we found that the inclusion of the four orientation angles does not alter the estimation of the binary masses and spins. 
The reason is that the orientation angles are rather uncorrelated with those parameters, appearing only in the GW 
amplitude and not in the phase. For the same reason, we found that spin-orbit and spin-spin couplings have little 
effect on the angular resolution, distance determination and hence on the redshift determination for massive black- hole 
binaries. For NS-IMBH binaries, these extrinsic parameters are determined rather poorly with or without spin effects. 
For massive black- hole binaries, Monte Carlo simulations show that LISA can provide reasonably accurate distance 
determinations out to z ~ 2 — 4 for black hole masses < 10 7 M Q . 

The cosmological reach of LISA will ultimately depend on design choices for the acceleration noise. The reason 
is that for massive binaries the GW signal sweeps through the low-frequency LISA band. By default we made a 
rather optimistic assumption that the LISA noise curve can be extrapolated down to a lower frequency of 10 -5 Hz. 
We then carried out an explorative survey to see how a higher (more conservative) low-frequency cutoff affects the 
determination of the binary parameters for high-mass configurations. We found that the cutoff will have a dramatic 
impact on our ability to i) locate black hole binaries of mass > 10 6 Mq in the sky, ii) measure their masses, iii) use 
them as standard cosmological candles, and iv) bound the mass of the graviton. Our results are compatible with 
similar investigations which have appeared in the literature |l8t Il9l. 

Our analysis was limited to non-precessing binaries. Vecchio [K)| has shown that for comparable high-mass BHs 
[e.g., (10 6 + 10 6 )Af Q ] modulational effects can decorrelate some of the binary parameters, allowing a better estimation 
of masses and distances with respect to the case when spins are aligned or antialigned with the (orbital) angular 
momentum. At this stage it is not clear if modulational effects can improve the accuracy in estimating binary 
parameters also for small mass-ratio binaries (e.g., a NS and an IMBH), and allow one to put more stringent bounds 
on alternative theories of gravity. Only a direct calculation can clarify this point, and we plan to tackle it in the near 
future. 

When describing stellar mass objects inspiralling into intermediate-mass black holes, we only considered circular 
orbits. This assumption is barely justified, especially for high mass ratios. In the future we plan to investigate how 
the results change when eccentricity is included. 

Finally, in this paper we have focused on statistical errors, implicitly assuming that the waveform is known with 
high enough accuracy to disregard systematic errors. However, massive binaries are likely to be detected with rather 
high SNR, on the order of 100. If this is the case, it might well be that spinning waveforms at 2PN order are not 
sufficiently accurate to permit one to neglect systematic over statistical errors. 
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APPENDIX A: ESSENTIAL FORMULAE FOR THE LISA RESPONSE TO NON-PRECESSING 

SPINNING BINARIES 

In this Appendix we write down the essential formulae of Ref. Q which we use in Sec. Ill Bl (We refer the reader 
to Ref. for further details and notation.) 

Unbarred quantities refer to the rotating LISA-b&sed coordinate system, while barred quantities refer to the fixed 
solar-system-based coordinate system. Assuming as in Ref. [2j that the noise is symmetric in each pair of LISA 
arms, we can reduce LISA to two independent Michelson interferometers with equilateral triangle shape. In this 
approximation the LISA beam-pattern functions for the two Michelson outputs are the same as for a single detector, 
except for the factor y/3/2 which already appears in Eq. (|2.23l) . and are given by 



h,(t>s,^s) = -(1 +cos 2 s )cos2(/)5Cos2?As 



cos 9s sin 2<ps sin 2ips ■> 



1 



,4>s,ips) = 2( 1+cos 0s)cos2^,ssin2V>s 
+ cos 9 S sin 2<ps cos 2ipg > 



and 



(Al) 



F*(9 s ,<l>s,rl>s) = F+(9 s ,<f> s --,il>s). 



(A2) 



In the above equations we have denoted by (9 s, 4>s) the source location and by ips the polarization angle defined as 

L • z — (L • n)(z • n) 



tan^s(t) = 



n • (L x z) 



(A3) 



L, z and — n being the unit vectors along the orbital angular momentum, the unit normal to LISA's plane and the 
GW direction of propagation, respectively. 

The waveform polarization and Doppler phases entering the GW signal l|2.23|) are (a — I, II): 



Vp,a(*) 



tan 



2tt/ 



2(L-n)F a *(f) 
(1 + (L • n)*)F+(t) 



= R sin 8s cos(<p(t) - (f> s ) , 

c 



(A4a) 
(A4b) 



with R — 1 AU and <j)(t) = c/)q + 2nt/T. Here T = 1 year is the orbital period of LISA, and (f>o is a constant that 
specifies the detector's location at time t — 0. In this paper we always assume that there is no precession, so L a 
points in a fixed direction (9l,4>l)- 

To express the angles (Qs,4>S,i>s) evaluated with respect to the rotating detector-based coordinate system as 
function of the angles (9s,$>£,9l,4>l) evaluated with respect to the fixed solar-system based coordinate system, we 
use the following relations 



1 - V3 
cos 9 s (t) = - cos 9 S — sin 9 S cos(<£(t) - <f> s ) 



h(t) = a + — + tan 



V3 cos 9s + sin9s cos(4>(t) — 4>s) 



(A5a) 
(A5b) 



2 sin 9s sm((f>(t) — <ps) 

where ao is a constant specifying the orientation of the arms at t = 0. Following Cutler 0, we take ao = and 
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4>q = 0, corresponding to a specific choice of the initial position and orientation of the detector. In addition, 

z • n = cos 9 S , (A6a) 
1 - V3 

L z = -cos 9 L — sin L cos (0(i) - 0l) , (A6b) 

L • n = cos ^lcos 9s + sin0£ sin ^5 cos (4>l — <f>s) , (A6c) 
n • (L x z) = - sin 9 L sin 9 S sin {<j) L - <j> s ) 

— -^-cos^(t) ^cos#l sin (9 s sin ^5 — cos (9~s sin^ sin0£ 
-f sin0 W (cosfc sin*, cos0 L - cos*, sin* s cos0 s ) . (A6d) 

APPENDIX B: SUBTLETIES IN THE INVERSION OF THE FISHER MATRIX 

In the paper we evaluate the Fisher matrix using both a Mathematica and a Fortran code. In the Fortran code 
we normally perform the numerical inversion of the Fisher matrix using the LU decomposition, which expresses the 
Fisher matrix as the product of a Lower-triangular and an Upper-triangular matrix (cf. Sec. 2.3 of |34|). To check 
the result we simply multiply the inverse by the original matrix. In this way we obtain a numerical "identity matrix" 
whose elements ij" um will be slightly different from the Kronecker symbol <5^. We can measure this deviation from 
the "true" identity matrix defining a small quantity 

ftasmaiy \llf m -6 %3 \ . (Bl) 

We found that extreme mass ratio inspirals are more likely to yield an ill-conditioned Fisher matrix. Therefore, as a 
rule of thumb, we consider the inversion successful if the parameter ej nv < 10~ 3 (for NS-BH binaries) and if ej nv < 1CP 4 
(for massive BH-BH binaries). In Mathematica we perform a similar check using the built-in matrix inversion routine. 

Matrix inversion generally becomes more difficult as the number of elements of the Fisher matrix increases. In 
particular, the LU decomposition and the Mathematica inversion routine fail when we consider alternative theories of 
gravity including both spin-orbit and spin-spin terms. To understand the reason for this failure we can use a principal 
component analysis |39f . also known as singular-value decomposition (see eg. Sec. 2.6 of [34j). We decompose the 
Fisher matrix F as 

F = UWV T , (B2) 

U and V being orthogonal matrices, and V T being the transpose of V. The matrix W is diagonal with positive or 
zero elements Wj (the singular values). The inverse of F is then given by 

F = VW _1 U T , (B3) 

where W _1 is a diagonal matrix whose elements are the reciprocals 1/wj of the singular values. Numerically speaking, 
a matrix is not invertible when the reciprocal of its condition number (defined as the ratio of the largest singular value 
to the smallest singular value) approaches the machine's floating-point precision. When we consider alternative 
theories including both spin-orbit and spin-spin terms, the matrix becomes non-invertible in this sense: our numerical 
experiments show that one of the singular values approaches zero. In principle, even in this case we can still obtain a 
"pseudo-inverse" : the matrix which is closest to the "real" inverse in a least-square sense 34] . To do this it suffices 
to replace 1/uij by zero whenever Wj is zero in Eq. I|B3I) . However, in this paper we decided not to quote results 
obtained in this way. We estimate the binary parameters only when spin-orbit and spin-spin terms ( Tables ITTT1 and IVj) . 
Brans-Dicke and spin-orbit terms ( Table 11*7]) or massive graviton and spin-orbit terms (Table |VT|) are included. 
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TABLE I: Number of GW inspiral cycles contributed by different PN orders for different NS-BH binaries. We assume 5 = 0.3 
and an observation time T b s = 1 yr. In the bottom section of the table, we normalize the number of cycles associated with 
the Brans-Dicke parameter to zu (first row) and to the Cassini bound wbd > i^Cassmi = 4 x 10 4 (second row). We also show 
the initial and final GW frequencies, assuming an upper cutoff of 1.0 Hz for the LISA noise curve. 



PN order 


(1.4 + 4OO)M 


(1.4 + 1000)M Q 


(1.4 + 5OOO)M 


(1.4 + 1O 4 )M 


MHz) 


4.601 x 10~ 2 


3.658 x 10~" 


2.446 x 10 _ - 


2.057 x 10~ 2 


/fln(Hz) 


1.000 


1.000 


0.8792 


0.4397 


Newtonian 


2,294,904 


1,828,036 


1,224,122 


1,025,711 


1PN 


35,366 


44,712 


67,309 


78,460 


Tail 


-18,064 


-29,081 


-66,278 


-89,793 


Spin orbit 


1,437/3 


2,314/3 


5,274/3 


7,145/3 


2PN 


422 


868 


3,016 


4,653 


Spin spin 


-139<r 


-288a 


-1,001a 


-l,545cr 


Brans-Dicke 


-3,560,569^ 


-1,793,782itj 


-536,954 m 


-319,126^7 


Brans-Dicke 


-89 OJCassini/l^BD 


-45 ^Cassini/^BD 


-13 OJCassmi/l^BD 


-8.0 WCassini/<^BD 
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TABLE II: Number of GW inspiral cycles contributed by different PN orders for high-mass BH binaries. We assume an 
observation time T b s = 1 yr. In the bottom section of the table, we normalize the number of cycles associated with the 
graviton- mass term to /3 9 (first row), and to the Compton wavelength A 9 , using Eq. 12.5bl (second row). We assume a 
luminosity distance Dl = 3 Gpc, Qm = 0.3 and J1a = 0.7. 



PN order 


(10' + 1O')M 


(10' + 10 b )M Q 


(10 b + 1O 6 )M 


(10 b + 10 4 )M Q 


(10 4 + 1O 4 )M 


fin (Hz) 


1.073 x 10" b 


2.361 x lO -5 


4.525 x lO -5 


4.199 x 10~ 4 


8.046 x 10 -4 


/fin (Hz) 


2.199 x 10" 4 


3.997 x 10~ 4 


2.199 x 10" 3 


3.997 x 10~ 2 


0.2199 


Newtonian 


535 


1,174 


2267 


21,058 


40,369 


1PN 


55 


115 


134 


677 


769 


Tail 


-48 


-127 


-92 


-450 


-308 


Spin orbit 


4/3 


10 p 


7/3 


36 P 


25/3 


2PN 


4 


8 


6 


18 


12 


Spin spin 


-1(7 


-2 a 


-Iff 


-5 a 


-3(7 


massive graviton 


-209/3 9 


-333/3 9 


-512/3 9 


-1967/3 9 


-2926/3 9 


massive graviton 


-1063 (l0 15 km/A 9 ) 2 


-478 (l0 15 km/A 9 ) 2 


-260 (l0 15 km/A 9 ) 2 


-28 (I0 15 km/A 9 ) 2 


-15 (I0 15 km/A 9 ) 2 



TABLE III: Errors and correlation coefficients for different NS-BH binaries at different PN orders in general relativity (no 
Brans-Dicke term) with and without spin-orbit and spin-spin terms. We consider one detector and set p = 10. 



PN order At c A<j> c 
(s) 


AM/M 

(%) 


Arj/n 
(%) 


A/3 


A(7 


c Mrl 


C M0 




c M ° 






(1.4 + 4OO)M 
1 1.59 3.61 


0.0000710 


0.0206 






-0.995 












1.5 


3.00 15.6 


0.000148 


0.149 






-0.999 












1.5 


4.07 28.1 


0.000478 


0.375 


0.00346 




-0.996 


0.951 


-0.918 








2 


3.33 17.6 


0.000220 


0.208 






-0.999 












2 


4.00 25.7 


0.000491 


0.393 


0.00260 




-0.996 


0.893 


-0.849 








2 


16.2 508 


0.00266 


3.85 


0.206 


1.91 


-0.996 


-0.981 


0.994 


-0.983 


0.995 


1.000 


(1.4 + 1000)M o 
1 1.86 3.41 


0.0000861 


0.0157 






-0.995 












1.5 


2.52 8.08 


0.0000233 


0.0369 






-0.922 












1.5 


4.45 25.2 


0.000560 


0.275 


0.00618 




-0.995 


0.999 


-0.991 








2 


2.55 7.89 


0.0000341 


0.0439 






-0.965 












2 


4.33 21.8 


0.000584 


0.297 


0.00560 




-0.996 


0.998 


-0.989 








2 


16.1 425 


0.00296 


2.65 


0.0705 


0.793 


-0.996 


-0.962 


0.982 


-0.980 


0.994 


0.997 


(1.4 + 5000)M o 
1 2.55 3.19 


0.000123 


0.0101 






-0.995 












1.5 


3.03 5.61 


0.0000386 


0.00870 






0.941 












1.5 


5.68 22.4 


0.000771 


0.168 


0.00673 




-0.995 


0.999 


-0.999 








2 


3.00 5.18 


0.0000360 


0.0101 






0.932 












2 


5.34 14.7 


0.000857 


0.202 


0.00700 




-0.996 


0.999 


-0.999 








2 


18.5 349 


0.00386 


1.52 


0.0137 


0.178 


-0.996 


0.951 


-0.920 


-0.975 


0.991 


-0.860 


(1.4 + 10 


4 )M 






















1 


3.33 3.43 


0.000154 


0.00914 






-0.995 












1.5 


3.96 5.82 


0.0000607 


0.00582 






0.965 












1.5 


8.60 26.7 


0.00103 


0.163 


0.00737 




-0.996 


0.998 


-0.999 








2 


3.91 5.26 


0.0000592 


0.00678 






0.964 












2 


7.72 13.0 


0.00120 


0.211 


0.00822 




-0.997 


0.999 


-0.999 








2 


34.6 487 


0.00576 


1.66 


0.0299 


0.133 


-0.997 


0.998 


-0.989 


-0.978 


0.992 


-0.961 
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TABLE IV: Errors and correlation coefficients in Brans-Dicke theory using 2PN templates, with and without the spin-orbit 
term. We consider one detector and set p = 10. In the first row we do not consider spin terms; in the second row we also include 
spin-orbit effects. When we include the spin-orbit term priors do not have an appreciable effect on parameter estimation. For 
each binary we also give the bound aj B D,unc that could be obtained (in principle) if all the binary parameters were known and 
not correlated with the BD term. 



At c A<p c AM/M At)/t) 
(s) (%) (%) 


ujbd A/3 




jW„ MB 






rod 


^ L.'-t -J - 4:UU J1V1Q , U>BlJ )Un c — 


43 057 645 












3.82 23.2 0.000243 0.293 


765,014 


-0.939 


0.421 


-0.705 






7.95 76.7 0.00657 2.50 


39,190 0.0508 


-0.997 


-0.997 0.999 


0.988 


-0.993 


-0.999 


(1.4+1OOO)M , WBD.unc = 


= 21,602,414 












3.79 16.7 0.000189 0.116 


211,389 


0.845 


-0.984 - 


-0.926 






7.99 58.4 0.00764 1.86 


21,257 0.0557 


-0.996 


-0.997 1.000 


0.987 


-0.998 


-0.995 


(1.4 + 5000)M Q , o;bd,u„c = 
4.60 12.5 0.000600 0.0342 


= 6, 388, 639 
50,925 


0.970 


-0.998 - 


-0.955 






8.79 23.4 0.0114 1.33 


6,486 0.0550 


-0.997 


-0.997 0.999 


0.988 


-1.000 


-0.992 


(1.4+1O 4 )M , WBD,unc = 


3, 768, 347 












6.59 13.8 0.000877 0.0253 


26,426 


0.979 


-0.998 - 


-0.963 






13.6 15.5 0.0178 1.61 


3,076 0.0706 


-0.998 


-0.998 0.999 


0.991 


-1.000 


-0.993 



TABLE V: Errors and correlation coefficients for various high mass BH binaries in general relativity using one detector, with 

and without spin-orbit and spin-spin terms. We set Dl = 3 Gpc and assume Ho = 72 km s _1 Mpc -1 . The effect of adding 
priors is practically negligible in all cases. 

~Af c A^ AM/M A^A? A/3 Act c Mv c Mli ?^ 

oo (%) m 

(10 Y + 10 Y )M Q , SNR=2063 

5.27 0.0108 0.00189 0.0401 - - 0.930 ----- 

7.26 0.0717 0.0224 5.56 0.315 - -0.996 0.996 -1.000 - 

79.6 1 12 0.0703 25.4 1.29 0.438 -0.994 0.997 -1.000 -0.948 0.976 -0.970 
(10 y + 10 S )M Q , SNR=1204 

5.11 0.0187 0.00104 0.0267 - - 0.934 ----- 

9.50 0.0589 0.0132 1.86 0.0891 - -0.996 0.997 -1.000 - 

83.4 2.18 0.0499 10.8 0.402 0.340 -0.995 0.999 -0.999 -0.964 0.985 -0.975 
(10° + 10 b )M Q , SNR=2143 

0.307 0.00551 0.000369 0.0157 - - 0.887 ----- 

0.496 0.00819 0.00303 1.02 0.0566 - -0.991 0.993 -1.000 - 

3.02 317 0.00872 4.60 0.213 0.140 -0.991 0.996 -0.999 -0.938 0.975 -0.964 
(10° + 10 b )M Q , SNR=2378 

0.214 0.00734 0.000202 0.00976 - - 0.871 ----- 

0.380 0.00773 0.00188 0.352 0.0162 - -0.991 0.994 -1.000 - 

2.46 521 0.00587 1.79 0.0560 0.0858 -0.991 0.999 -0.995 -0.948 0.980 -0.957 
(10 s + 10 S )M Q , SNR=1710 

0.0521 0.00539 0.000114 0.0113 - - 0.799 ----- 

0.0746 0.00614 0.000766 0.366 0.0196 - -0.985 0.989 -1.000 - 

0.300 0.201 0.00210 1.70 0.0685 0.0892 -0.987 0.995 -0.997 -0.931 0.977 -0.958 

(W b + 10 4 )M Q , SNR=601 



0.122 0.0165 0.0000810 0.0111 - - 0.732 ----- 

0.190 0.0301 0.000767 0.225 0.00934 - -0.989 0.994 -0.999 - 

0.643 0.663 0.00219 1.08 0.0213 0.101 -0.988 0.995 -0.970 -0.937 0.978 -0.899 

(10 4 + 10 4 )M Q , SNR=252 

0.222 0.0419 0.0000822 0.0308 - - 0.846 ----- 

0.327 0.103 0.000622 0.612 0.0275 - -0.984 0.991 -0.999 - 

0.758 1.39 0.00200 3.39 0.0786 0.480 -0.990 0.998 -0.984 -0.951 0.984 -0.937 
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TABLE VI: Errors and correlation coefficients for various high mass BH binaries including the massive graviton and spin-orbit 
terms. We use one detector, set D L — 3 Gpc and assume H = 72 km s _1 Mpc -1 . We show the massive graviton bound 
obtained assuming (Qm = 0.3, SIa = 0.7). The effect of adding priors is practically negligible in all cases. For each binary we 
also give the bound A 9 . unc that could be obtained (in principle) if all the binary parameters were known and not correlated 
with the massive graviton term. 



~Kf c Af c AM/M A^A? K A/3 c Mr] c Mli « c r,p " c Mp c 73 ^" 

(s) (%) (%) (10 15 km) 

(TO 7 + W')M Q , SNR=2063, A 9 , unc /(10 i!) km) = 880 

5.27 0.0108 0.00189 0.0401 - - 0.930 ----- 

14.1 0.0448 0.0155 0.534 69.4 - -0.981 -0.993 0.997 - 
79.6 1.12 0.0703 49.2 22 2 3.08 -0.978 0.948 -0.994 0.975 -1.000 0.995 
(10 Y + 1O B )M , SNR=1204, A 9 , un c/(10 ib km) = 527 

5.11 0.0187 0.00104 0.0267 - - 0.934 ----- 

13.8 0.0749 0.0104 0.352 39.5 - -0.985 -0.995 0.997 - 
83.4 2.18 0.0499 25.6 9.57 1.52 -0.981 0.964 -0.997 0.978 -1.000 0.998 
(10 b + 1O B )M , SNR=2143, A 9 , unc /(10 ib km) = 351 

0.307 0.00551 0.000369 0.0157 - - 0.887 ----- 

0.675 0.0175 0.00244 0.146 46.3 - -0.968 -0.988 0.994 - 
3.02 0.317 0.00872 12.2 12.4 0.790 -0.963 0.938 -0.997 0.960 -1.000 0.997 
(10° + 1O S )M , SNR=2378, A 9 , unc /(10 i:> km) = 215 
0.214 0.00734 0.000202 0.00976 - - 0.871 - 

0.481 0.0238 0.00161 0.0865 27.4 - -0.973 -0.992 0.994 - 
2.46 521 0.00587 5.52 6.02 0.337 -0.966 0.948 -0.998 0.962 -1.000 0.999 
(10 b + 10 b )M o , SNR=1710, A 9 , unc /(10 ib km) = 139 

0.0521 0.00539 0.000114 0.0113 - 0.799 - 

0.0843 0.0131 0.000678 0.0713 23.2 - -0.952 -0.986 0.987 - 

0.300 0.201 0.00210 6.58 4 92 0.436 -0.950 0.931 -0.998 0.947 -1.000 0.999 

(10 s + 1O 4 )M , SNR=601, A^ncAlO 1 '" km) = 77.4 

0.122 0.0165 0.0000810 0.0111 - - 0.732 

0.203 0.0458 0.000712 0.0804 10.8 - -0.973 -0.994 0.990 - 

0.643 0.663 0.00219 5.48 1 76 0.352 -0.950 0.937 -0.999 0.946 -1.000 1.000 

(10 4 + 1O 4 )M , SNR=252, A 9 . unc /(10 lb km) = 48.5 

0.222 0.0419 0.0000822 0.0308 - - 0.846 

0.333 0.120 0.000597 0.217 5.80 - -0.964 -0.990 0.990 - 

0.758 1.39 0.00200 29.7 0.670 2.06 -0.957 0.951 -1.000 0.955 -1.000 1.000 
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TABLE VII: Average errors for a Monte Carlo simulation of 10 4 NS-BH binaries randomly located and oriented on the sky 
with mass of (1.4 + 1000)Mq. We first consider general relativistic waveforms (GR) and add spin-orbit (SO) and spin-spin 
(SS) couplings, fixing the single-detector SNR pi — 10. Then we do the same including also a Brans-Dicke (BD) term. In each 
case, the first line refers to the errors obtained using only the first detector; the second line gives the errors obtained using both 
detectors; the third line gives pattern-averaged results from the relevant entries of Tables If III and ll VI 



Case 


AD L /D L 


AM/M 

(%) 


Arj/rj 

(%) 


Afl s 
(10~ 5 str) ( 


59s 
arcmin) 


WBD 


A/3 


Act 


GR 


0.782 
0.376 


0.0000633 
0.0000378 
0.0000341 


0.103 
0.0617 
0.0439 


15.3 
5.95 


42.5 
26.5 


- 


- 


- 


GR+SO 


0.797 
0.374 


0.00178 
0.00111 
0.000584 


0.816 
0.507 
0.297 


35.3 
13.8 


64.6 
40.4 




0.0179 
0.0111 
0.00560 




GR+SO+SS 


0.900 
0.420 


0.00869 
0.00562 
0.00296 


6.96 
4.51 
2.65 


60.2 
23.9 


84.3 
53.1 




0.143 
0.0930 
0.0705 


1.82 
1.18 
0.793 


BD 


0.764 
0.359 


0.000789 
0.000488 
0.000189 


0.283 
0.174 
0.116 


40.8 
15.8 


69.4 
43.2 


62,561 
96,719 
211,389 






BD+SO 


0.898 
0.413 


0.0225 
0.0143 
0.00764 


4.87 
3.10 
1.86 


86.8 
34.1 


101.3 
63.5 


7,209 
10,799 
21,257 


0.157 
0.100 
0.0557 





TABLE VIII: Average SNRs and errors for a Monte Carlo simulation of 10 4 BH-BH binaries randomly located and oriented in 
the sky with mass (10 6 + 1O 6 )M0. We fix Dl = 3 Gpc and, where indicated, include a massive-graviton (MG) term, spin-orbit 
(SO) and spin-spin (SS) couplings. We also illustrate the deterioration in parameter estimation when we assume that LISA is 
blind below some cutoff frequency /i ow . The default value for /i ow is 10~ 5 Hz. In each case, the first (second) line refers to the 
errors of a (10 6 + 1O 6 )M0 binary using one (two) detectors. We assume a cosmological model with Qm = 0.3 and Qa = 0.7. 



Case 


SNR AD L /D L 


AM/M 


Av/V 


AQ S 


(Az/z) 


(Az/z) bcst 


A 9 


A/3 


Act 








(%) 


(96) 


(10 -5 str) 






(10 15 km) 






GR 


1861 


0.0458 


0.000614 


0.0304 


59.9 


0.101 


0.0380 










2693 


0.0106 


0.000349 


0.0155 


9.76 


0.0873 


0.00880 








GR+cutoff 


1823 


0.341 


0.00966 


0.109 


5894 


0.306 


0.283 








/low = 10" 4 Hz 


2640 


0.0164 


0.00121 


0.0289 


48.9 


0.0881 


0.0136 








GR+SO 


1863 


0.0492 


0.00558 


1.93 


64.4 


0.103 


0.0408 




0.107 






2696 


0.0107 


0.00295 


1.01 


10.2 


0.0873 


0.00891 




0.0558 




GR+SO+SS 


1862 


0.0504 


0.0152 


8.05 


67.1 


0.104 


0.0418 




0.374 


0.248 




2695 


0.0109 


0.00852 


4.52 


10.4 


0.0873 


0.00902 




0.209 


0.139 


MG 


1861 


0.0486 


0.00447 


0.273 


64.1 


0.103 


0.0403 


37.4 








2693 


0.0107 


0.00237 


0.145 


10.2 


0.0873 


0.00889 


49.5 






MG+cutoff 


1787 


0.442 


0.0443 


1.18 


13290 


0.386 


0.367 


16.2 






/low = 10" 4 Hz 


2592 


0.0159 


0.00921 


0.362 


51.2 


0.0878 


0.0132 


29.5 






MG+cutoff 


1859 


0.0511 


0.00522 


0.301 


71.8 


0.104 


0.0424 


35.4 






/ low = 5 x 10" 5 


Hz 2691 


0.0108 


0.00274 


0.158 


10.6 


0.0873 


0.00895 


46.7 






MG+SO 


1861 


0.0495 


0.0152 


21.5 


67.2 


0.103 


0.0411 


10.6 


1.39 






2693 


0.0108 


0.00852 


12.1 


10.4 


0.0873 


0.00896 


13.3 


0.780 
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FIG. 1: Analytic approximation to the LISA root noise spectral density curve used in this paper and in Ref. fl3| (dashed line) 
and the curve produced using the LISA Sensitivity Curve Generator 14] (solid line). The SCG curve has been multiplied by a 
factor of \/3/20 to obtain an effective non-sky averaged noise spectral density (see Sec. Ill CJI . The SCG noise curve does not 
include the extragalactic white dwarf confusion noise while the analytical approximation curve does. 




FIG. 2: Left: luminosity distances Dl of NS-BH binaries observed with SNR=10 as a function of the black hole mass. We 
assume that the NS mass Mns = 1.4M©. Right: SNR for equal mass BH-BH binaries at Dl = 3 Gpc as a function of the total 
mass. Solid lines refer to the LISA noise curve 112.311 used in this paper; dashed lines refer to the same noise curve without 
including the white-dwarf confusion noise. The "bump" in the noise curve due to white-dwarf confusion noise is responsible for 
the dip in the SNR for MBH binaries of masses ~ 10 6 M Q . 
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FIG. 3: Monte Carlo simulation of 10 4 binaries with observed total mass (1.4 + 1O 3 )M0 in general relativity, with single-detector 
SNR=10, £7a = 0.7, Qm = 0.3. Top: probability distributions of the angular resolution Afls in steradians for one detector 
(left) and two detectors (right). In each figure, from left to right, the histograms refer to no spins, SO included, and SO and 
SS included. Bottom: probability distributions of ADl/Dl for one detector (solid) and two detectors (dashed); ADl/Dl is 
essentially unaffected by the inclusion of spins, so we only show histograms without the spin terms. 
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FIG. 4: Monte Carlo simulation of 10 4 binaries with observed total mass (1.4 + 1O 3 )M0 in general relativity, with single-detector 
SNR=10, = 0.7, Qm = 0.3. Top four panels: probability distribution of the errors on the chirp mass AA4/M, the reduced 
mass Afi/fi, the SO parameter A/3 and the SS parameter Act. Bottom panel: bound on o>bd when a Brans-Dicke term is 
included. Solid (dashed) lines refer to one (two) detector(s). 
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FIG. 5: Monte Carlo simulation of 10 4 binaries with total mass (10 6 + 10 6 )Mq in general relativity, with Dl = 3 Gpc, Qa = 0.7, 
fljlf = 0.3, with no spins. Panels show probability distributions of the SNR, the distance determination error ADl/Dl, the 
redshift errors Az/z and (Az/z)b CBt and the angular resolution Afls in steradians. Solid (dashed) lines refer to one (two) 
detector(s). 
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FIG. 6: Monte Carlo simulation of 10 4 binaries with observed total mass (10 6 + 1O 6 )M0 in general relativity, with Dl = 3 Gpc, 
= 0.7, Qm = 0.3. Top four panels: probability distributions of the error on the chirp mass AM/M, the reduced mass 
A (j, /(j,, the SO parameter A/3, the SS parameter Act. Bottom panel: bound on the graviton Compton wavelength A 9 (in km), 
when a massive graviton term is included. Solid (dashed) lines refer to one (two) detector(s). 
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FIG. 7: Errors on different parameters as a function of the redshift z for nonspinning BH-BH binaries with source mass 
(10 6 + 10 6 ) Mq (solid lines) and (10 7 + 10 7 ) M Q (dashed lines). We fix fi A = 0.7, Qm = 0.3. Errors are obtained by averaging 
results of Monte Carlo simulations of 10 4 binaries and using two detectors. Left: general relativistic non-spinning binaries 
(GR); right: general relativistic binaries including SO and SS terms (GR+SO+SS). Circles refer to chirp mass M, squares to 
reduced mass ft, triangles to angular resolution £ls in steradians, diamonds to luminosity distance Dl- 




FIG. 8: Errors on different parameters as a function of the total binary mass M (in solar masses) for equal-mass nonspinning 
BH-BH binaries. We fix Dl = 3 Gpc, S7a = 0.7, Q.m = 0.3. Errors are obtained by averaging results of Monte Carlo simulations 
of 10 4 binaries and using two detectors. Left: general relativity (GR); right: massive graviton theories (MG). Solid lines assume 
that the LISA noise can be extrapolated down to /i ow = 10 -5 Hz; dashed lines assume /i ow = 10 -4 Hz. Circles refer to chirp 
mass M, squares to reduced mass fj,, triangles to angular resolution fig in steradians, diamonds to luminosity distance Dl- 
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FIG. 9: Bound on the graviton Compton wavelength as a function of the total binary mass M (in solar masses) for nonspinning 
equal-mass BH-BH binaries. We fix Dl = 3 Gpc, S1a = 0.7, SIai = 0.3 and use both detectors. The solid line assumes that the 
LISA noise can be extrapolated down to /i ow = 10~ 5 Hz; the dashed line assumes /i ow = 10~ 4 Hz. 



